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ABSTRACT 


A scheme is introduced which makes it feasible to make completely 
self-consistent Brueckner-Hartree-Fock (BHF) and renormalized BHF calcu- 
lations for spherical, closed-shell and axially-symmetric deformed nuclei. 
The usual requirement or orbital (HF) self-consistency has been imposed, 
as well as self-consistency in the starting energies and occupation 
probabilities. Previously only approximate forms have been used for the 
Pauli operator. This approximation is removed here and a method for mak- 
ing the necessary Pauli corrections to the reaction matrix during the 
approach to self-consistency is presented. A discussion of the symme- 
tries which reduce the problem to one of manageable proportions is in- 
cluded. 
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SELF-CONSISTENT PAULI CORRECTIONS IN BRUECKNER-HARTREE-FOCK CALCULATIONS 


by Richard C. Braley and William F. Ford 
Lewis Research Center 


SUMMARY 

The Brueckner-Hartree-Fock method is extended to include self- 
consistency in the reaction matrix as well as the usual self- 
consistencies. Corrections to the G-matrix are made using an approximate 
iterative procedure rather than the exact matrix inversion. Considera- 
tion of only spherical and axially symmetric deformed nuclei makes it 
possible to reduce the computational problem to one of manageable propor- 
tions by making use of the symmetries of the nuclear wave functions. 


INTRODUCTION 

Nuclear many-body theory has now progressed to the point where it 
becomes possible to make a quantitative connection between the (realistic) 
nucleon-nucleon interaction and properties of finite nuclei. These theo- 
retical developments (referred to generally as Brueckner theories) are 
excellently summarized in several recent review articles (refs. 1 to 3). 
Applications of the theory to spherical (closed-shell) nuclei have been 
numerous (refs. 4 and 5) , and recently the theory and its application 
have been extended to deformed nuclei (ref. 6). Significantly, all of 
these Brueckner-type calculations have yielded saturation at radii which 
are too small. 

On the other hand, ordinary Hartree-Fock (HF) calculations but with 
density-dependent or velocity-dependent effective interactions have been 
very successful in predicting binding energies and shapes of finite nu- 
clei (refs. 7 and 8). Despite the considerable degree of phenomenology 
inherent in such investigations , it is apparent that their underlying 
assumptions are good approximations to the realities of the nuclear 
medium, and that there is much to be learned from them. In particular, 
those who employ the more complicated Brueckner theories must wonder 
whether some neglected aspect of these theories could provide the equiva- 
lent of density or velocity dependence. 

So far, Brueckner-type calculations have been carried out in a 
largely but not completely self-consistent manner. The usual requirement 
of "orbital" or Hartree-Fock self-consistency has been imposed, as well 
as self-consistency in the starting energies and (in renormalized BHF 
calculations) in the occupation probabilities; but only approximate forms 
have been used for the Pauli operator, which forbids scattering into in- 
termediate states which are normally occupied. Since the correct form 
for the Pauli operator is not determined until the calculation has been 
completed, all current methods have the drawback that their output is in- 
consistent with their input. While this drawback is recognized and de- 
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plored, the enormous computational difficulties involved have so far pre- 
vented any attempts to treat the Pauli operator exactly. It is certainly 
conceivable, however, that current approximations make the resulting ef- 
fective interaction (reaction matrix) less sensitive to details of the 
nuclear environment than would otherwise be the case, and that the suc- 
cess of the phenomenological calculations lies in their simulation of 
this sensitivity by means of density or velocity dependence. 

If complete Pauli corrections are sought, it will be necessary to 
make corrections to the reaction matrix during the approach to self- 
consistency. A method for doing this in the general case of a deformed 
nucleus will be presented here, together with a discussion of those sym- 
metries which reduce the problem to one of manageable proportions. 


METHOD 

All modern theories of nuclear matter and finite nuclei make use, in 
one form or another, of the reaction matrix G(to) which satisfies 

G(m) = V + V Q G(u>) (1) 

This equation describes the scattering of two nucleons which are moving 
in some average field, with V representing the interaction between them 
and a) the "starting energy." The unperturbed motion of the nucleons 
(i.e., in the absence of V) is described by the two-body Hamiltonian h, 
which is usually taken to be a sum of one-body Hamiltonians 

h = h^ + 


in which T represents the kinetic energy operator and U the average 
field felt by the nucleon. 

The reaction matrix is used in place of V in many-body theories, 
the argument being that in this way multiple scattering between a given 
pair of nucleons is treated exactly. For this reason it is important to 
exclude, from the description of the multiple scattering, any effects 
which have already been included in the description of the unperturbed 
nuclear state. This is accomplished by means of the projection operator 
Q, which prohibits either of the nucleons from scattering into any of the 
single-particle orbitals which make up the unperturbed ground state. Be- 
cause this restriction is similar to that arising from the indistinguish- 
ability of the nucleons, the operator Q is often called the Pauli oper- 
ator . 


Evidently the Pauli operator, and hence the reaction matrix, cannot 
be determined unless the occupied single-particle orbitals are known. 
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But in all the theories , the single-particle orbitals are eigenstates of 
a Hamiltonian which depends in some way on the reaction matrix. Thus yet 
another type of self-consistency is required, which has been avoided in 
calculations so far by the use instead of an approximate reaction matrix 
G 0 (uj) obtained with an approximate Pauli operator Q 0 . (It should be re- 
emphasized that this approximation is only a computational expedient, 
since in principle the correct Q could be obtained by solving for the 
reaction matrix and the nuclear orbitals simultaneously.) 

Actually the approximation goes farther than that, because the two- 
body Hamiltonian h also depends on the single-particle orbitals. 

Usually it is replaced by a "reference" Hamiltonian, h Q , and the corre- 
sponding G 0 (u>) is termed the "reference" reaction matrix. 

In any event, G(co) and G 0 (<o) satisfy a relation (ref. 1) of the 

form 

g(“) - V“> + - rvTr) GW < 2 > 

derived by eliminating the two-body interaction V. This equation has 
been used in some instances to obtain corrections to the reference reac- 
tion matrix, but not in a fully self-consistent way. It is possible to 
achieve complete self-consistency, but the corrections must be made more 
than once as the calculation unfolds. 


Since the determination of nuclear structure is basically a nonlinear 
problem, all methods of solution which do not approximate the nonlinearity 
must proceed iteratively. Thus at any stage, say the n^, there is a 
single-particle Hamiltonian H( n ) , whose eigenfunctions may be labeled 

In the iterative process these eigenfunctions are used to form the 
next Hamiltonian H( n+ 1), and the calculation continues until the eigen- 
functions and eigenvalues do not change from iteration to iteration. The 
various methods in current use (HF, BHF, RBHF) differ from one another in 
practice only in the way in which the new Hamiltonian is generated from 
the previous eigenfunctions. 


To date all such calculations have employed a fixed reaction matrix 
determined in some approximate manner at the outset (ref. 9). But a 
slight modification of (2) provides a means of correcting the reaction 
matrix at each iteration. For, once the eigenfunctions have been 

found, the Pauli operator QW corresponding to them may be constructed, 
and the corrected reaction matrix G^ n ) obtained by solving 


G<n> 



+ G*^ 
o 



(n) 

T<s> 


co 



(n) 


(3) 
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G (n) = G (n-1) + G (n-1) 

and used to form the new Hamiltonian as before. When convergence 

is near, it is legitimate to replace G( n ) by G( n- 1) on the right-hand 
side of (4) and then the reaction matrix changes at each iteration by the 
(matrix) amount 

AG = X 1 - X Q (5) 

NEW 

X 1 (m) = G(u) q — G(o») (5a) 

oi - n 




,(n-l) 


M h (n) 
u) - n ■ 


m - h 


(n-1) 


,(n) 


(4) 


.OLD 


X Q ((o) = G(o>) 


ui - h 


OLD 


G(u>) 


(5b) 


SYMMETRY CONSIDERATIONS 

All current nuclear structure calculations require that the reaction 
matrix be expressed in terms of a suitable single-particle basis. Gener- 
ally the basis functions are spherical, i.e., they may be specified by 
the quantum numbers (n, Z, j, m, t}, Z representing the orbital angular 
momentum, j the total angular momentum, m the z-component of total 
angular momentum, and x the z-component of isospin, with n standing 
for any other quantum numbers needed to distinguish the state (usually 
the radial quantum number). In the coordinate representation, a partic- 
ular basis function labeled by the index "a" will be denoted by 

E E ^ n a l a j aVa) (6) 

of course standing for the space and spin coordinates of the nucleon. 

The reaction matrix is a two-body operator, and so its matrix ele- 
ments will be of the form ^ab | G(u) | cdy. Because only antisymmetrized 
matrix elements are ever needed, 

^b|G(aj)|cd^ = -^ab ( G (to) jdc^ = -^ba |G(io) [ cd^ (7) 

and it suffices to consider pairs A = {ab} and B = {cd} such that 
a < b and c < d. In most calculations the reaction matrix is real, and 
it is always Hermitian (for fixed w) , so that 

<a|g(w)|b) = ^b|g(u,)|a) 

Thus only matrix elements with A <_ B need to be calculated. The net 


( 8 ) 
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result of all these simplifications is to reduce the number of matrix 
elements required by essentially a factor of 8. 

Further reductions can only come about from symmetry effects, such 
as invariance with respect to spatial and isospin rotations and time- 
reversal invariance. Unfortunately, although V always possesses such 
symmetries, G(w) may not because of its dependence on the nuclear orbi- 
tals through h and Q. For instance, although all stationary states of 
nuclei are also states of sharp angular momentum, the so-called intrinsic 
state of a permanently deformed nucleus - for which the basic structure 
calculations are performed in one version of the current theory - is 
characterized by having a rather broad spread of angular momenta forming 
a basis for the rotational "band." 

In like manner, the isospin and time-reversal symmetries possessed 
by stationary nuclear states are generally not found in the intrinsic 
state. The exception occurs for doubly-even nuclei with N = Z, whose 
orbitals may be filled in such a way as to preserve the symmetries. Fur- 
thermore, these nuclei are quite often axially symmetric, so that the 
z-component of total angular momentum is conserved. The enormous reduc- 
tion in calculation time which results from all this makes such nuclei 
extremely desirable for theoretical studies. 

Accordingly, it is convenient to divide the pair states A = {ab} 
into groups with the same parity (tt^ = ir a + tt^) , the same z-component 
of angular momentum (M A = m a + nq,) , and the same total isospin 
(T a = 0 or 1, t a = x a + Tfc). In this way the reaction matrix is made 
block diagonal, each block being characterized by the appropriate con- 
served quantum numbers. 

When isospin is conserved a further simplification results, since 
the T A = 1 matrix elements are independent of t a (a consequence of 
the Wigner-Eckart theorem). Thus only two of the four submatrices need 
to be computed. 

A similar result occurs when time-reversal invariance holds, but for 
different reasons. To begin with, the pair states of the basis may al- 
ways be arranged so that for each state |a) there is a time-reversed 
counterpart |S^), the two satisfying 

I*) ■ 31 A>. |A> - a|A> (9) 

where is the time-reversal operator. The matrix element of an arbi- 
trary operator Q between |T) and |"B^ may then be written (ref. 10) 

(a|£2|b)= (^A^I^b) (10a) 

= ^J7A | 27 B) (10b) 

= (aI^qJb)* 


(10c) 
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<^A | « | b)> = <A | Q" _1 Q S7 | b)>* = (A|ft|B>* (lOd) 

where ft is the time-reversed form of ft. (The step from (10b) to (10c) 
follows from the antiunitary nature of the time-reversal operator.) When 
the operator ft is time-reversal invariant, 

(a|a|b)* = <A | ft | b)* (11) 

and when the phases have been chosen so that the matrix element is real, 

^A|ft|B)= ^A | ft | b) (12) 

Note that ft = ft only for nuclei whose intrinsic state is time reversal 
invariant . 

One immediate consequence of (12) is that no negative-M submatrices 
need be computed, since each will be equal to the corresponding positive- 
M submatrix. For M = 0 the situation is more involved, and it is con- 
venient to further subdivide this set so that one group of states is 
mapped into the other by time reversal. There are many ways of doing 
this; one simple method is to place pair states with m a < m^ in the 
group [A] , and those with m a > m^ in_the group [A] . By the foregoing 
argument matrix elements of the type ^A | ft | A */ need not be computed, and 
when the pair states are from different groups, 

^A|ft|A') = ^A|ft|A’) (13a) 

= ^A ' | ft | a) (13b) 

(Eq. (13a) follows from (12), and (13b) from (8).) By symmetry, then, 
itv is only necessary to compute 

^A|ft|A'^, A <_ A* , and (14a) 

^A|ft|A'), A_<A' (14b) 


when M = 0. 


DETAILS OF THE CALCULATION 


Both terms in (5) are similar in form, so that what must be calcu- 
lated are matrix elements of the general type 



(15) 


Now in all current versions of nuclear, many-body theory, the eigenfunc- 
tions 4> ^ may be used to construct (in a manner to be specified shortly) 
orbital pair states | Ay) which diagonalize both Q and h. Hence the 
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needed matrix elements may be obtained from 


f 



( 16 ) 


where l_j indicates the legacy of Q, namely that only the unoccupied 
orbitals are to be included in the summation. 


This expression may be simplified in several ways. First, because 
of the antisymmetric nature of the matrix elements, the summation may be 
restricted to those orbital pairs A = {Xy} with X < y and the result 
then multiplied by 2. Second, if the orbital pairs are labeled by the 
same general quantum numbers that label the basis pairs, the summation 
will be further restricted to those pairs with the same quantum numbers 
as {ab} and {cd}. Equation (15) may thus be. written 


Y <“> A 



G(m)^G(m)^ 
co - e. 


(17) 


where Y^y indicates the further restrictions and an obvious shorthand 
for the matrix elements has been introduced. 

The economy of expression displayed in (17) is achieved only through 
the use of mixed-representation matrix elements of G(co). These are re- 
lated to the basis-representation matrix elements by 


G (^ 2 £ G( »)K (18) 

B 

where is the projection of |A^ onto the basis pair-state |b^: 

|a )-C p b | b> (l9a) 

B 

Pg = (b|A> (19b) 

(The factor of 2 appears as before because of the use of ordered pairs.) 
Although calculation of the mixed-representation matrix elements is a 
rather lengthy intermediate step which must be performed for both terms 
in (5) , overall time is greatly reduced owing to the appearance of a 
single, rather than a double, pair-state summation in (17). 
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To proceed further it is necessary to introduce an explicit notation 
for the isospin. Let x x (12) denote singlet (T = 0) and triplet (T = 1) 
isospin functions for two nucleons 
with a definite isospin is then of 

■ ^(<?l |a K2l b ) 

= ^x .2 |abT^x T (12) 


| abT ) = L- { | a > 1 b> + (-1) T | b ) ( a > } (20c) 

V2 

A similar decomposition of the orbital pair state | X^ can be made, lead- 
ing to the expression 



or symbolically 


labeled 1 and 2. A pair state |ab) 
the form 


+ (-D T 4 1 |b)<^ 2 |a)jx x (12) 


(20a) 


(20b) 


p aS = C A + <-!> C b C a 


where is the projection of | X^ onto the basis state | a^ : 


( 21 ) 


i^> - X) c - U > 

a 

C a ' 4 |X ) 

In this notation two-body matrix elements are written 

(abT|G(to) |cdl) = (ab|G T (u>) | cd) = G T (co)^ 


(22a) 


(22b) 


(23) 


and have the exchange character (-1)^. (It should be noted that the pro- 
jection of the time-reversed orbit | x) onto the time-reversed basis state 
]a^ is also equal to C^, provided this coefficient is real.) 

The only remaining problem, aside from the technical questions of 
efficient sequencing and indexing, is the energy-dependence or co- 
dependence of the reaction matrix. It is clear from its definition that 
G(co) is a meromorphic function with poles at the eigenvalues but 

in practice these poles, which form the intermediate-state spectrum, are 
expected to be relatively far removed from the range of co-values needed 
for boundstate investigations. Such an assumption must be verified in, 
each instance, of course, but when true it permits an accurate approxima- 
tion of the form G(oo) 3 %(to) with 
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G n (w) 


N 


D 

n=0 


G P (o>) 

n n 


(24) 


where P n (w) is a polynomial of order n. Equation (24) is Inconvenient 
for Pauli corrections when the basis space is very large, and so it is 
rearranged to read 


v u) 


N 


Yj G N (t \> )/? v (aj) 
v=0 


(25) 


where each ^(u>) is a Lagrange-interpolation polynomial of order N 
characterized by 

-fit* ) = 6 (26) 

v p pv 

With this formulation the Pauli corrections may be performed directly on 
the matrix coefficients %(ui v ). The roots are quite arbitrary, 

but round-off error is minimized if they are confined to the m-range of 
interest and spread rather uniformly throughout it. The coefficients 
g n( u v) may be obtained directly from G(m) , if desired, or a least- 
Bquares procedure used to obtain (24) which is then evaluated at the 
various values 


CONCLUDING REMARKS 

The problem of making Pauli and spectral corrections to the reaction 
matrix for spherical and axially-gymmetric deformed nuclei has been dis- 
cussed. This additional self -consistency in BHF and RBHF calculations is 
usually avoided due to the complexity of the numerical computations. A 
significant reduction in complexity of the problem results when an itera- 
tive method is used to obtain the new G-matrix (instead of matrix inver- 
sion) , and the symmetries of the intrinsic state are used to reduce the 
number of matrix elements which must be computed. The completely self- 
consistent RBHF problem is then reduced to one of manageable proportions. 
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